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FOREWORD 


Since  the  number  of  artificial  earth  satellites  requiring  operational  flight 
control  is  continually  increasing,  more  emphasis  has  been  placed  upon  the  develop¬ 
ment  of  rapid  and  accurate  ephemeris  prediction  techniques.  This  study  has  been 
conducted  under  the  auspices  of  the  Defense  Mapping  Agency  in  order  to  develop  a 
formalism  that  can  be  adapted  to  aid  in  the  satisfaction  of  operational  control 
requirements  associated  with  low-altitude  earth  satellites.  This  report  was 
reviewed  and  approved  by  Dr.  R.  J.  Anderle  and  Mr.  R.  W.  Rill. 

Released  by: 

o  £  /^LfLA^ 

THOMAS  A.  CLARE,  Head 

Strategic  Systems  Department 


DTIO 


*^ELECTE 

AUG  1  4  885 


B 


Accession  For 

1  NTTS  GRA&I 

*  DTIC  TAB 

□ 

Unannounced 

□ 

1  Justification.  — - 

T itt 

Dy 

Distribution/ 

Availability 

Codes 

[Avail  end/or 


Diat  i  Special 


NSWC  TR  84-353 


•u*» 


Section 

1 

2 

3 


APPENDIXES 

A 

B 


CONTENTS 

Page 

INTRODUCTION .  1 

EQUATIONS  OF  MOTION  IN  NONSINGULAR  VARIABLES .  2 

DRAG  DECAY  RATES  FOR  NONSINGULAR  ELEMENTS  IN  AN 

OBLATE  DIURNAL  ATMOSPHERE  .  3 

GEOPOTENTIAL  PERTURBATIONS  IN  NONSINGULAR  ELEMENTS .  10 

REFERENCES .  17 

POSITION  PARTIAL  DERIVATIVES  WITH  RESPECT  TO  ORBITAL  ELEMENTS  .  18 

A  DERIVATION  OF  THE  INTEGRAL  EXPRESSIONS  FOR  Aa»  A£,  An  .  .  .  .  20 
DISTRIBUTION .  (1) 


NSWC  TR  84-353 


1 .  INTRODUCTION 


There  has  in  recent  years  been  a  renewed  interest  in  satellite  operational 
control  using  rapid  and  accurate  prediction  techniques.  Of  especial  interest  are 
satellites  in  low  earth  orbits  that  experience  atmospheric  drag  disturbances. 

Early  studies  concerning  the  effects  of  atmospheric  drag  perturbations  upon 

1  2 

artificial  satellites  of  small  eccentricity  were  performed  by  Cook  and  Kin*  'ele,  * 

3  4  5 

Zee,  Otterman  and  Lichtenfeld,  and  Lane.  More  recent  analyses  have  been 

6  7  8 

developed  by  Liu  and  Alford,  and  Santora.  *  The  latter  author's  works  represent 

one  of  the  more  realistic  approaches  to  the  determination  of  drag  perturbations  on 

9 

a  decaying  satellite.  There  the  oblateness  and  diurnal  characteristics  of  the 
atmosphere  have  been  combined  with  the  effects  induced  by  the  oblate  earth  mass 
gravitational  perturbations  to  provide  a  "unified”  analytic  theory  describing  both 
the  drag  and  gravitational  influences  exerted  upon  a  low  altitude  earth  satellite. 

Santora' s  theory  is  unified  in  the  sense  that  the  oblate  diurnal  atmospheric 
drag  perturbations  have  been  combined  with  those  produced  by  the  terrestrial 
gravity  field.  It  is  not  a  completely  unified  theory  for  low  altitude  satellite 
motion  because  a  significantly  different  development  must  be  used  to  treat  orbits 
with  0.  S  e  <  .01  than  is  used  for  orbits  with  .01  <  e  s  .20.  A  completely 
unified  theory  is  presented  here  that  not  only  combines  the  oblate  gravitational 
effects  with  those  produced  by  an  oblate  diurnal  atmosphere,  but  also  eliminates 
the  eccentricity  dichotomy  that  exists  in  Santora 's  approach. 

The  following  sections  of  this  report  discuss  in  detail  the  development  of  this 
unified  theory.  The  nonsingular  element  set  and  the  associated  Lagrange  planetary 
equations  used  in  the  development  of  the  theory;  a  detailed  derivation  of  the  drag 
decay  rates;  and  an  overview  of  the  geopotential  perturbations  are  presented  in  the 
following  sections.  Appendixes  are  also  provided  that  contain  mathematical 
relationships  used  in  the  development  of  the  decay  rate  expressions. 
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2.  EQUATIONS  OF  MOTION  IN  NONSINGULAR  VARIABLES 

The  eccentricity  dichotomy  that  exists  in  the  Santora  theory  can  be  removed  by 
redeveloping  the  theory  in  terms  of  the  following  element  set: 


a  *  semimajor  axis 
A  =  M  +  a)  +  ft 
5  ■  e  cos  <L  (&  =  (d  +  0) 
n  =  e  sin  w 
P  =  sin  (i/2)  cos  ft 
Q  *  sin  (i/2)  sin  ft 


(2.1) 


where  a,  e,  i,  w,  ft,  and  M  are  the  usual  Keplerian  elements.  This  element  set  is 
nonsingular  for  V  /  rr  and  e  <  1. 

The  Lagrange  planetary  equations  for  the  nonsingular  element  set  defined  in 
Equations  2.1  are: 

2 

a  =  —  R, 
na  X 

J  •  "  -  -h  Ra +  tjl 2«rc  +  "V  +  rV(PRP  +  <V 

2na  2na  y 


2?}~~  {RX  -  ^2  Rn  -  TV  "<PRP  +  «V 
na  (1  +  y)  na  2na  y 


i/i  ~  ?<PRp +  QV 

na  (1  +  y)  na  2na  y 


P(nRc  -  5Rn) 


(2.2) 
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where  n  is  the  mean  motion  and 


Y  *  >1  ~  e  . 


(2.3) 
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For  close-earth  artificial  satellites  perturbed  by  the  earth's  oblateness  and 
atmospheric  resistance,  the  partials  of  the  disturbing  function  are  of  the  form 


R 

a 


(2.4) 


where  a  represents  any  of  the  elements  of  Equations  2.1,  and  the 
and  G  designate  the  perturbing  terms  due  to  atmospheric  drag  and 
oblateness,  respectively.  To  be  more  explicit 


superscripts  D 
the  earth’s 


(2.5) 


and 


G  =  3R 
a  3a 


(2.6) 


where  f  is  the  nonconservative  drag  force,  r  is  the  satellite  radius  vector,  and 
R  is  the  gravitational  disturbing  function. 


3.  DRAG  DECAY  RATES  FOR  NONSINGULAR  ELEMENTS  IN 
AN  OBLATE  DIURNAL  ATMOSPHERE 

Analytic  expressions  for  drag  decay  rates  in  an  oblate  diurnal  atmosphere 
using  nonsingular  variables  are  developed  in  this  section.  Since  the  principal 
drag  effects  on  a  close-earth  satellite  are  perturbations  in  a  and  e,  it  is 
necessary  only  to  consider  expressions  for  the  drag  induced  changes  in  a,  £, 
and  n. 

An  analytic  form  for  the  atmospheric  density  can  be  obtained  by  combining 
the  oblate  atmosphere  used  by  Cook,  King-Hele,  and  Walker^  with  that  of  an 
atmosphere  with  diurnal  variation  discussed  by  Cook  and  King-Hele.^  The  resulting 
form  for  the  density  is  given  by 

p  *  Pq(1  +  F  cos  <}>)exp[-B(r  -  o)]  (3.1) 

where  Pq  is  an  average  atmospheric  density  defined  by 


3 


(3.2) 


P0  2^pmax  +  Pmin^ 


F  is  a  density  amplitude  factor  given  by 


F  = 


P  -  P  . 

max  mm 

P  ,  P  . 

max  +  mm 


(3.3) 


B  =  H 


0  =  r 


-1 


,  .2  .2 

1  -  e  sin  x  sin  u 


.  .2  .2 

1  -  e  sin  x  sin  w 


(3.4) 

(3.5) 


and 


cos  <p  =  A 


cos  E  -  e 


1  -  e  cos  E 


+  B 


(1  - 


2n  2  • 

e  )  sin 


1  -  e  cos  E 


(3.6) 


In  the  last  four  equations  p  and  p  .  are  the  maximum  daytime  and  minimum 

max  min 

nighttime  densities,  respectively,  at  the  osculating  perigee  altitude  h  ;  H  is 
the  density  scale  height;  e  is  the  earth's  ellipticity;  r^  is  the  osculating 
perigeal  radius;  E  is  the  eccentric  anomaly;  u  =  6  +  w,  where  0  is  the  true 
anomaly ;  and 

A  =  sin  6  sin  i  sin  u  +  cos  6  <  cos  (£  -  a„)  cos  m 

D  D  \  D 


-  cos  i  sm 


(fi  -  0^)  sin  wj 


(3.7) 


and 


sin  sin  x  cos  u>  -  cos  6  <  cos  (fi  -  aD)  sin 
B  D  |  B 

+  cos  i  sin  (JJ  -  ag)  cos  wj 


(3.8) 


The  declination  and  right  ascension  of  the  center  of  Ihe  diurnal  bulge  are 

designated  by  <5  and  aD,  respectively,  in  the  last  two  equations. 

B  B 


Equation  3.5  may  be  expanded  to  first  order  in  to  give 
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1  2 

o  =  rpU  +  ~2  e  sin  i  (cos  2u  -  cos  2w)]  (3 

Similarly,  Equation  3.6  may  be  expanded  to  first  order  in  e  to  give 

2 

cos  $  =  A(cos  E  -  e  4-  e  cos  E)  +  B(sin  E  +  e  sin  E  cos  E)  (3 

Substituting  Equations  3.9  and  3.1J  into  Equation  3.1  and  using  the  relation 

r  =  a(l  -  e  cos  E)  (3 

allows  the  following  first-order  expression  to  be  written  for  the  atmospheric 
density: 

2 

p  =  Pq[1  +  FA(cos  E  -  e  +  e  cos  E)  +  FB(sin  E 
+  e  sin  E  cos  E)]exp|  -6ae(i  -  cos  E) 

+  c  cos  2u  -  c  cos  2u;|  (3 

where 

1  2 

c  =  -=:  eBr  sin  i  (3 

2  P 

This  quantity  may  be  treated  as  a  small  parameter  of  the  same  order  of  magnitude 
as  the  eccentricity  so  that  the  following  expansion  may  be  used: 


/  l  12  2 

exp  <  c  cos  2u  |  =  1  +  c  cos  2u  +  —  c  cos  2u 


The  changes  in  the  elements  a,  and  n  over  one  orbital  revolution  due 
to  atmospheric  drag  deceleration  are  given  by  (see  Appendix  B) 


e  cos  E)~ 
e  cos  E) 
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/CDS\  ,  f  ( \  + 

-  I  -jjj-  I  bay  I  pi  -j— I 


AC  = 


1 

+  e 

cos 

1 

-  e 

cos 

e] 

/2 


(, 


CC'i  Ui  CO  d  t. 


-  sin  &  sin  E)  dE 


(3. 


and 


_/VV  r27l  4  e  cos  cV 
\  m  /  aj  c  yi  -  e  cos  e y 


An  = 


/  2 


(y  sin  i  cos 


+  cos  Si  sin  E)  dE 


(3.1 


where 


(¥). 


is  the  inverse  ballistic  coefficient  for  the  satellite,  and 


6  =  (1  -  — 2-  Aw  cos  t )  ^ 
v  e 
P 


(3.1 


In  the  last  expression  is  the  speed  of  the  satellite  at  perigee  and  A  is 

the  ratio  of  atmospheric  angular  rotation  rate  to  that  of  the  earth,  designated  by 

W  . 

e 


Using  the  relations 


cos  0  = 


cos  E  -  e 
1  -  e  cos  E 


(3.1' 


and 


1/2 


sin  0  = 


(1  -  e  ) _ sin  E 


1  -  e  cos  E 


(3 . 2< 


to  eliminate  0  to  first  order  in  e;  applying  Equations  3.12  and  3.14;  and 
expanding  the  resulting  integrands  of  Equations  3.15  -  3.17  yields: 


Aa 


■  ■(¥)••• 


a  exp[-8ae  -  c  cos  2 id]  j  (I  4-  FA  cos  E 

0 


{<■ 


+  eFA  cos  E  -  eFA  +  FB  sin  E  +  eFB  sin  E  cos  Ej • <1 

2 


+  c  cos  2  (a.  +  E)  -  2ce  sin  2  (^  -t-  E)  +  — 
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12  2  )  ( 
+  —  c  cos  4(w  +  E)  -  c  e  sin  4(u>  +  E)  sin  E>*<1 

+  2e  cos  e|  exp  [Bae  cos  E]  dE 


=  -  dp^ayexp [-Bae  -  c  cos  2u)]^  |sAMe|  •  |sAMe|  •  |i 

+  e  cos  e|  •  |y  cos  w  cos  E  -  sin  ui  sin  e|  expfBae  cos  E]  dE  (3. 


An  =  - 


(c  s\  /*  ^ 

- J  dp^ayexp^Bae  -  c  cos  2 u,]j  |sAMe|  ■  |sAMe|  •  |l 
+  e  cos  e|-|y  sl-n  w  cos  E  +  cos  w  sin  e|  exp  [Bae  cos  E]  dE  (3. 


When  the  integrands  of  the  last  three  equations  are  multiplied,  the  results 
contain  trigonometric  terms  that  are  expressible  as  functions  of  cos  (nE) ,  n  =  0, 
1,  2,  7.  This  permits  them  to  be  written  in  terms  of  the  integral  represen¬ 

tation  of  the  Bessel  function  of  the  first  kind  and  imaginary  argument  defined  by 


I  (Bae) 
n 


*/o 


cos  (nE)exp(Bae  cos  E)  dE 


The  changes  in  a,  i,  and  n  over  one  orbital  revolution  due  to  atmospheric  drag 
deceleration  then  become 


Aa  =  -  2tt 


/CDS\  2  ( 

tt  I  1  dp^a  exp[-Bae  -  c  cos  2w]  <  1^  +  2el^  +  FA[I^ 

+  |  e(IQ  +  3I2)  ]  +  c  |l2  +  2el3  +  ~  FAHj  +  I3 

+  j  e(IQ  +  2I2  +  51^)]  |  cos  2oo  -  ~  cFB{Ii  ~ 

+  |  e(I0  +  4I2  -  5I4)  }  sin  2o>  +  ^  c2  j  IQ  +  2el1 


+  FA[I1  +  \  e(I0  +  3I2)]|  +  |  c2|l4  -  e(I3  -  3I5> 


If 

K5 
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+  J  FAtlj  +  I5  -  |  e(I9  -  21  ^  -  /'f  :  ]icos  4a» 

- 1  -  h  - 1  e«2  ^  e\  T  sLn  4 


/  v  \  i 

— y  6pQaexp[~8ae  -  c  cos  2u.l  <  j  i  I  +  —  i\\e ) I ^ 
+  |  (e  +  FA)(IQ  +  I  )  +  4-  FAeL3](l  +  V)  cos  £ 


(3.25) 


-  -  FB[Iq  -  I2  +  e(I1  -  I3)](l  +  >-■)  air.  .. 

+  {f  [I1  +  X3  +  FAe ( 1 3  +  I5>  -  T  e(lQ  -  2I2  -  3I4) 

+  Y  FA(Iq  +  2I2  +  1^)  ]  cos  £  +  |-  FBcf  1^  -  2I2  +  1^ 

+  2e(I1  -  2I3  +  I5)]  sin  ioj  cos  2c  -  -  c  |fb[Iq 

-  I4  +  2e(I1  -  I5)]  cos  u  +  [e(I0  -  4I2  +  31^) 

+  2(I3  -  I1)  -  FA(I0  -  X4  +  2e(I3  -  I,))]  sin  0  J  sin  2w 
c2  ( 

+  16  { [2(I3  +  I5)  +  e(5I6  +  2I4  -  3I2)  +  FA(16  +  2I4 


+  I2  +  e(3l?  +  3I5  -  I3  -  1^)]  cos  u,  +  FB [I  -  2I4 

+  I2  +  e(3l^  -  71^  +  5I3  -  i. ) ]  sin  \  |  cos  4c 
c2  / 

+  (  FB[I6  -  I2  +  e(3I?  -  I5  -  313  +  1  j)  J  cos  S 


+  [ 2 ( I 3  -  l5)  -  e(5I6  -  8I4  +  3I2)  -  FA(Lfc  -  I, 
+  e(3Iy  -  3I5  -  I3  +  Ij ) ) ]  sin  £  }  sin  4  j 


(3.26) 
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res' 

An  =  — 2tt  j6pQaexp[-8ae  -  c  cos  2o>]^-|  FB(1  +  -  1 2 


'{i 


+  e(I1  -  I3)]  cos  S3  +  [ <1  +  \  eFA)I1  +  y(e  +  FA)(I0 
+  I2)  +  2  eFAI3Hl  +  Sin  S3  -  FB[I0  -  2I2  +  I4 
+  2e(I1  -  2I3  +  I5)]  cos  u  -  ■|[I1  +  I  +  eFA(I3  +  I5> 

+  e(“  J  Io  +  X2  +  I  V  +  FA(I  Xo  +  I2  +  I  V1  sin  “}  cos  203 

+  f  {[e(I0  -  4I2  +  3I4)  -  2(IX  -  I3>  -  FA(IQ  -  I4  +  2e(I3 

-  I^))]  cos  S3  -  FB[Iq  -  I4  -  2e(I^  -  1^)]  sin  sin  2uo 
c2  ( 

‘T6\FB[I6  ~  2I4  +  X2  +  e(3I7  "  7I5  +  5I3  "  ll)] 


i 


vVi 

m 


COS  0) 


-  [2(I3  +  I5)  +  e(5I6  +  2I4  -  3I2)  +  FA(I6  +  2I4  +  I2 
+  e(3l?  +  3I5  “  I3  “  Ij)) ]  sin  S3  J  cos  4w  -  j[2(I3  -  I5> 


-  e(5l6  -  8I4  +  3I2)  -  FA(I6  -  I£  +  e(3I?  -  3I5  -  I3 


+  IL) ) ]  cos  S3  -  FB[Ig  -  I2  +  e(3ly  -  I$  -  3I3 


+  1^)]  sin  S3 


|  sin  4w| 


(3.27) 


It  should  be  noted  that  y  does  not  appear  in  the  last  two  equations  since  it 
has  been  replaced  with  its  first-order  expansion  equivalent  of  unity.  The 
changes  in  eccentricity  and  orbital  period  per  revolution  can  be  found  from 


Ae  -  (AC2  +  An2) 


1/2 


(3.28) 


and 
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At  =  4  -  Aa  (3.29) 

l  a 

respectively,  where  t  is  the  orbital  period. 

4.  GEOPOTENTIAL  PERTURBATIONS  IN  NONSINGULAR  ELEMENTS 

In  order  to  evaluate  the  drag-induced  changes  in  a,  and  n  given  by 
Equations  3.25  through  3.27  on  a  rev-by-rev  basis,  it  is  necessary  to  propagate 
the  gravitationally  perturbed  six-member  nonsingular  element  set  to  the  time  of 
interest.  This  subsection  is  devoted  to  a  development  of  the  gravitational 
disturbing  function  R  of  Equation  2.6. 

As  will  be  discussed  in  subsequent  studies,  computations  performed 
using  this  development  will  be  done  in  terms  of  the  mean  state.  Thus,  it  is 
convenient,  yet  completely  general,  to  divide  the  gravitational  disturbing 
function  into  short-periodic  and  secular-long-periodic  parts: 


R  -  RSLp  4-  Rsp  (4.1) 

where  the  subscripts  "SLP"  and  "SP"  denote  the  secular-long-periodic  and 
short-periodic  parts,  respectively,  and  are  given  by 

l  l 

RSLP  =  Z  Z  Z  Z  R£mpq 

t> 2  m=0  p=0  q=2p-£ 

and 

l  l 

vEE  I  E  r*.pi  (4-3) 

t> 2  m=0  p=0  qj*2p-£ 

The  (^mpq)1"^1  terms  in  the  disturbing  function  summations  are^ 


1 


ftl 
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(£) 


,<A*«  cos  6£mpq 


+  hm  3i"  Vpq>  +  ^tmpq'V,  sln  Vopq 


‘  B*n  cos  6*npq> 


where  ag  is  the  mean  equatorial  radius  of  the  earth,  p  is  the  earth's 


gravitational  constant 


Clm’ 

£  -  m 

even 

_S£m 

,  £  -  m 

odd 

S£m’ 

£  -  m 

even 

C£m’ 

i 

a 

odd 

0£mpq  =  U  '  2P  +  <1>X  -  m0 


Here  0  is  the  Greenwich  sidereal  time. 


The  J^p  and  K£p£J  functions  in  Equation  4.4  are  the  inclination  and 
eccentricity  functions  given  by 


Jj?_n  (c)  =  (-Dk 
Imp  2£p,( 


g  +  m) !  y  l  nj/2l-2p\/  2p  \  2£-a-2j 

.«-,)•  til  v j  /v-r 


o  j+(a-|a|)/2 
(1  -  c  ) 
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"  |q|+k  k  f  r  /  „  \ 

%q(v,  -  NA,  +  S  £  £  m-  te.) 

k=0  r=0  t=0  ' 


(1  +  Y)r+t*k(l  -  Y)k  ,  (for  q  >  0) 


%q<Y) 


£  tb  L&r) 

k=0  r=0  t=0L  '  ' 


(2p:«)^y 


(1  +  y)r+t"k(l  -  Y)k  »  (for  q  <  0) 


where 


c  =  cos  (i/2) 

k  =  integral  part  of  [^p] 
j  ^  =  max(0,  -  a) 
j2  =  min(2£  -  2p,  l  -  m) 
a  =  m  -  l  +  2p 


The  -Rfmoa  and  Ufa  functions  in  Equation  4.4  are  given  by 


k  U2 


n=0  u=u . 


(4.9) 


(4.10) 
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I&npq  =  Z 


n=0 


(-l)n+u+1S 


u=u^' 


U(  iai  \ 

u\u  JUn+l-u) 


|q|  -u 


^Up  |a|  -2n-l+Up2n-u+l 


(4.12) 


where 


k  -  [-M  ±  Ml,  k'  .  jjal  *  i°i  -  ij 

=  max(0,  2n  -  |a|),  =  min(2n,  |q|) 

=  max(0,  2n  +  1  -  |a|),  =  min(2n  +  1,  |q|) 


6^  =  1,  if  q,  a  are  both  positive  or  negative 


Su  -  (-1)  ,  if  q  or  a  is  negative 


The  partial  derivatives  of  the  inclination  and  eccentricity  functions  that 
are  required  for  use  in  the  Lagrange  planetary  equations  are 


and 


3J 
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3P 
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3  J  »  3  J  n 
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3Q 
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3K-  .  3K„ 

ipq  =  _  I  ipq 

3?  y  3y 


3K»  3K» 

ipq  =  _  a  ipq 

3r\  Y  3Y 


(4.13) 


where 
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The  partial  derivatives  of  the 


Rand  I 


»  (q  <  0) 


functions  are  as  follows 
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"N 
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Impq 


3£ 


=  |q| 
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£mpq' 


/  ’  3n 
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(4.17) 


3 

3P 


^mP3  =  |a|  ¥  •  '^Pq  =  JR 

11  ^m'pq’  3Q  JLX^' 


J 


where 


q  -  1  ,  (q  >  0) 
q  +  1  ,  (q  <  0) 


and 


/ 

m 


m  -  1  ,  (a  >  0) 
m  +  1  ,  (a  <  0) 


The 

standard 


nonsingular  element  set  of  Equations  2.1  can  be  decomposed  into  the 
Keplerian  element  set  through  application  of  the  following  relations: 
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and 


2  2 
e  =  U  +  n  ) 


1/2 


-12  2 

i  =  2  sin  [(P  +  Q  )  ] 


-1 


a  =  tan  (Q/P) 


cj  =  tan  1  (ti/£)  -  Q 


M  =  X  -  (to  +  n) 


(4.18) 


In  closing  it  should  be  mentioned  that  the  osculating  perigee  radius  r^  is 
required  to  compute  the  decay  rates  of  Equations  3.25  -  3.27.  This  quantity 
is  given  by 


rp  rSLP  +  rSP 


(4.19) 


evaluated  at  0  =0.  The  rgLp  portion  is  computed  using  the  a  and  e  obtained 
from  the  Lagrange  planetary  equations  when  the  gravitational  disturbing  function 
R  is  represented  by  RgLp  only.  The  short  periodic  variation  of  the  satellite 
radius  is 


(4.20) 


where  p  is  the  orbital  semi-latus  rectum. 
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APPENDIX  A 

POSITION  PARTIAL  DERIVATIVES  WITH  RESPECT  TO  ORBITAL  ELEMENTS 


Relationships  that  are  required  for  the  transformation  of  R^  into  a  more 
usable  form  are  collected  in  this  appendix.  Let  x  j  =  1,  2,  3,  represent  the 
inertial  Cartesian  coordinates  of  the  satellite  position  vector  r.  Then 


3a  a 


3x  3x 

-  =  —  [—  - - e(sin  oi  cos  0  +  cos  oi  sin  ft  cos  i)] 

3X  y  r  3u 


3x2  ,  3x2 

-  =  —  [—  - - e(sin  ui  sin  ft  -  cos  oi  cos  ft  cos  i)] 

3X  y  r  3u 


a  1  3X3 

— —  =  —  [—  — —  +  e  cos  oi  sin  i] 
3X  y  r  3u 


by 


-  =  -r[sin  u  cos  ft  +  cos  u  sin  ft  cos  i] 

3u 
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=  r  cos  u  sin  \ 
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APPENDIX  B 

A  DERIVATION  OF  THE  INTEGRAL  EXPRESSIONS  FOR  Aa,  AC,  An 

In  order  to  obtain  integral  expressions  for  Aa,  AC,  and  An,  it  is  necessary 
to  develop  the  Gauss  planetary  equations  from  the  Lagrange  planetary  equations  of 
section  2  for  the  case  when 


R  =  R  =  F 
a  a 


This  is  accomplished  by  first  resolving  the  disturbing  force  F  into  three 
components  U,  V,  and  W,  where  W  is  perpendicular  to  the  plane  of  the  orbit 
(positive  toward  the  north  pole);  V  is  in  the  plane  of  the  orbit  at  right 
angles  to  the  radius  vector  (positive  in  the  direction  of  motion);  and  U  is 
along  the  radius  vector  (positive  in  the  positive  r  direction) .  The  required 
transformation  is: 


cos  u  cos  ft 

-  sin  u  sin  ft  cos  i 


=  cos  u  sin  ft 

+  sin  u  cos  ft  cos  i 


F  sin  u  sin  v 

z 


-  sin  u  cos  ft 


sin  ft  sin  i  U 


-  cos  u  sin  ft  cos  i 


-  sin  u  sin  ft 


-  cos  ft  sin  i  V  (B. 


+  cos  u  cos  ft  COSl 


cos  u  sin  i 


Applying  Equations  B.2  along  with  the  relationships  of  Appendix  A  to 
Equation  B.l  gives  the  following  results  for  R^: 

R»  -  <f>“ 


R°  =  (—  sin  6)U  +  (£-*■) V 


V 
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the  U,  V,  W  representation  to  the  T,  N,  W  representation,  where  T  is  the  component 
along  the  velocity  vector;  N  is  normal  to  T  (positive  toward  the  interior  of  the 
ellipse);  and  W  is  normal  to  the  orbital  plane  as  before.  The  required  transfor¬ 


mation  is: 


(1  +  e  cos  8) 


e  sin  0 


e  sin  6 


-  (1  +  e  cos  0) 


where 


a  =  -y  1  +  e  +2e  cos  0 

.  * 

Applying  this  transformation  to  the  U,  V,  W  representation  for  a, 
and  n  and  retaining  only  the  T  component  yields: 


2T 

a  =  —  o 

ny 


£  =  — —  [ 2  cos  (5  +  0)  +  e(sin  (5  +  8)  sin  8 
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n  =  —^-[2  sin  (5  +  0)  +  e(sin  (5  +  8)  cos  0 


-  cos  (5+0)  sin  0  +  sin  5)] 


The  retarding  force  due  to  atmospheric  drag  resistance  is 


w. 


'  ti  > 
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where  5  is  given  by  Equation  3.18  and  p  is  the  atmospheric  density.  Substi¬ 
tuting  Equation  B.18  into  Equations  B. 15  -  B.17;  using  Equations  3.19  -  3.20; 
and  applying  the  relations 


v2  -(v \(  1  +  e  cos  E 

\  a  J  \  1  -  e  cos  E 


n  > 

\jL 

1  -  e  cos  E^ 

f  dE 

allows  the  following  to  be  written 
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-  e  cos  E  I 


The  changes  in  these  elements  over  one  orbital  revolution  are  obtained  by 
integrating  the  last  three  equations  over  2tt  radians  of  eccentric  anomaly  to 
give  the  final  results 
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